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ABSTRACT 

Obtaining a confidence interval for a parameter X of an exponential 
distribution is a frequent occurrence in life testing problems. Often- 
times the test plan used is one in which all the observations are censored 
at the same time point t . Several approximate confidence interval pro- 
cedures are available in the statistical literature; however, to the 
knowledge of the author, the performance characteristics of the various 
approximations used in these procedures have not been established 
analytically. The purpose of this paper is to report the results of an 
empirical study of the performance of four of these procedures with 
respect to the expected length of the interval, the variance of the 
interval length, and the coverage probability. 
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I. INTRODUCTION 



In life testing applications, it is frequently desired to obtain a 
confidence interval for the parameter \ of an exponential distribution. 

In case a test plan is used for which all the observations are censored 
at the same time point t Q , several approximate confidence interval pro- 
cedures are available in the statistical literature. To the knowledge 
of the author, the goodness of the various approximations used in these 
procedures, and hence the goodness of the procedures themselves, have 
not been established analytically. The purpose of this paper is to 
report the results of an empirical study of the performances of four of 
these procedures with respect to the expected length of the interval, 
the variance of the interval length and the coverage probability. 

The general setting of the problem is as follows: suppose the random 

variables T^, > • * • > are independent and identically distributed with 

common exponential distributions 

f i^i^ = e t i > 0; 1 = 1 » 2> •*’’ n - 

The random variables that will actually be observed are X^, X^ y . . . , , 

where 



X. = 
1 



T. 

l 


if 


T. 

l 


< 


t 

o 


t 

o 


if 


T. 

l 


> 


t ; 1 = i, 2 > . • • > n 



These random variables are sometimes referred to as ’’censored” exponential 

variables, and t^ is called the ’’censoring point”. Define the random 

variables Y 1 , Y 0 , . . . Y by 
12’ n 



Y. = 

l 



1 if X. < t 
l o 



if X. > t ; i = 1 , 2 , ...,n. 
l—o 
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It is easily seen that the random variables , . . . , Y^ are independent 

and have Bernoulli distributions with parameter p = 1-e ^o. In Section II, 
four confidence interval procedures for X based on the and Y^ are 
discussed . 
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II. CONFIDENCE INTERVAL PROCEDURES 



A. PROCEDURE I 

This procedure is obtained as a special case of a solution to a more 

general problem that has been derived by Halperin [1], The random 

variable Y = has a binomial distribution with unknown parameter p. 

Suppose y is the observed value of the random variable Y, then the 100 

(1-CX )$ confidence bounds p and p may be obtained by solving the 

u ti 

equat ions 

y 

P [Y < y] = 2 (f) P L (l-p)"" 1 = a/2 

i=o 

and 



P [Y > y] = 2 Qp 1 (l-p) nU = a/2 

i=y 

for p. Due to the discreteness of the random variable Y, these equations 

do not yield an exact 100(l-a)$ confidence interval. The confidence 

coefficient is at least 100(l-a)$. The confidence interval for X is 

easily computed by inverting the relation p = 1-e ^o to get 

\ L = -ln(l-p L )/t o and = -ln(l-p u )/t o . 

A useful tool in solving for p and p is the Incomplete Beta 

L U 

function , 



F(x;a,B) - ^11 f t a (1 . t) S dt ml _ £ (^ + 1A i (1 . x) a + B + l-i 

* P * J i=o 

o 

The upper confidence bound, p^, is obtained by solving 

F^; y.n-y-1) = l-a/2 ; 
and the lower bound, p^, is obtained from 

f (p l ; y.n-y-1) = a/2. 
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The simplicity of the computation involved in this method is very 
apparent; however, it should be noted that the observed values of the 
random variables X^, X 2 , . . . , X^ have been disregarded. Intuitively 
it would appear that this procedure will not produce as good results 
as one which incorporates all the information of the experiment. How- 
ever, due to its simplicity, this procedure merits consideration. 

B . PROCEDURE 2 

For X t < < 1, Y = £Y. is nearly a Poisson random variable with 
parameter X2X^. Rubenstein [2] shows that^ = 

i 

nearly unbiased estimator for X, where the second factor is used to 

correct for the bias of the otherwise minimum variance estimator 

Y/£X . Wilks [3] states that for a Poisson random variable P with 

parameter p, the following relationship is approximately true: 

p-P = Zs/p where Z is the standard normal random variable. Rubenstein 

applies this reasoning to the random variable Y along with his modified 

estimator^ to obtain the relationship X - ^ = 2aJx . This equation is 

solved for X n and X by using the appropriate standard normal percentage 

point corresponding to a 100(l-a)$ confidence level. The resulting 

solutions for X and X are 

r ^ 9 9 A 9 

X L = 1_2X + Z C - (4tz C + Z C ) J/2 

and 

1/2 

X U = [ 2 ^ + 7,2 C + ^ z2q + Z^) ] /2 where C = ( a i) 1/2 ‘ 

C . PROCEDURE 3 

This procedure is a modification of one which has been suggested by 
Birnbaum [4] for application to uncensored life testing problems. Using 
the terminology commonly used in the literature of life testing to 



[l + j-] is a 
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facilitate the explanation, imagine that the random variables X^, , 

. .., X^ are observed sequentially. That is, a randomly selected item 
is put on test and is replaced with a similar item at failure or after 
a period of time t Q has elapsed, whichever occurs first. Then the 
arrival process of failures is a Poisson process. Suppose the experi- 
ment is continued until a total of k failures have been obtained. It 

is well known that the individual inter-arrival times of failures are 

t h 

exponentially distributed and the time to the k failure has a gamma 

distribution. Consider a test plan in which the experiment is stopped 

after a random amount of time (as would be the case for example, if n 

items were to be tested). The number K of failures is then a random 

variable. However, it would appear that, given K=k, the distribution 

of the time, W^, until k failures have arrived can be approximated by a 

gamma distribution. More precisely, suppose the observed value of the 

k k 

random variable K = E Y. is k, and let W, = E X.. Then the con- 

. , 1 k . , 1 

1=1 i=l 

ditional distribution of W^, given that K=k, is approximately the gamma 
distribution 

( . X k k-1 -Xw 
f (w) = p-^y v e ; w > o. 

It follows that V = 2XW can be approximated by a Chi-square variable 

2 2 

with 2k degrees of freedom. Thus, if x l-Ct/2 and X Ot/2 are the upper 
and lower a/2 percentage points of the Chi-square distribution with 2k 



degrees of freedom, then 



o 



X 2 a/2 y 2 l-oJZ' 



2W. 

k 



2W, 



2 "\ 

- ) constitutes an approximate 

K 

lOO(l-a)# confidence interval for X. 

Birnbaum also suggests an estimator for X which is merely the average 
of the upper and lower confidence bounds. This estimator has also been 
computed and tabulated. 
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D. PROCEDURE 4 



In Procedure 3 the distribution of the random variable W , the waiting 

K. 

t h 

time until the k failure, is approximated by a gamma distribution. 

Since the test is terminated after n items are tested, the maximum value 

that can attain is nt^. Consider the conditional probability 

P [W, < w and W. < nt ] P [W. < w] 
k — k — o _ k — 

P [W. < nt ] ~ P [W t < nt ] * 

k — o k — o 

This yields a truncated gamma distribution having density function 



P [W. < w|W- < nt ] 
k — 1 k — o 



f (w) 

pY . 7 o < w < nt 

F (nt ) — — o 



f T (w)= 



1 . 



w < o, w > nt 



Where f(w) is the density function of the gamma distribution and 

nt 

r ° 

F(nt^) = | f(w) dw. Intuitively, this new function would seem to 



o 

approximate the probability distribution of W^ in a censored test. An 
obvious drawback of this method is the difficulty of computation since 
both X and F(nt^) which depends on X are unknown. Thus, given tables 
of the Chi-square distribution function one is forced to use an itera- 
tive method for obtaining the desired X such that 

P [x 2 < 2\w, ] P [X 2 < 2Xw ] 

S — = a/2 5 — = l-a/2 

P [x < 2\nt ] P [x < 2Xnt ] 

o o 

to obtain \ and \ respectively. 

L U 

The estimator for X described in Procedure 3 is applicable to this 
procedure and is also tabulated. 
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III. COMPARISON OF PROCEDURES 



A Monte Carlo study has been made to compare Procedures 1, 2, and 3 
described above. One thousand samples of size n (n = 30, 40, 50) from 
an exponential distribution with parameter X (X = .1, .2, .8, 3, 5, 10) 
have been generated. For each sample, confidence intervals for X have 
been obtained by using each procedure for various censoring points t^ 
and confidence coefficients (1-a) (a = .05, .01). The tabulated 
quantities are the average length of the confidence interval, the sample 
variance of these lengths, and the empirical coverage probability of the 
intervals (i.e., the proportion of intervals which actually covered X ) . 

An abbreviated Monte Carlo study has been made to compare Procedures 
3 and 4. The lengthy computation of Procedure 4 required that the number 
of repetitions be reduced to 100 and that fewer combinations of X, t^, n, 
and a be used . 
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IV. CONCLUSIONS 



The results of this study differentiate each procedure as to its 
merits and shortcomings with respect to certain life testing situations. 
Each method is discussed below in order to define the situations in which 
it could be used. The following discussion includes only comparisons for 
Procedures 1, 2, and 3. Procedures 3 and 4 are compared separately and 
to a lesser degree due to the differences in the Monte Carlo studies made. 

Procedure 1 performed as expected; it generally gave less accurate 
results than the other procedures with respect to all three quantities-- 
expected interval length, sample variance of this length, and empirical 
coverage probability. However, since the empirical coverage probability 
of this procedure tends to be conservative and since the computation 
required for this procedure is minimal, Procedure 1 would be favored in 
cases wher. a quick but dependable confidence interval is needed or when 
a rough estimate is needed for use in more sophisticated procedures. 
Disregarding the observed value of the random variable is the main 

reason for the conservative results. 

The overall performance of Procedure 2 seems to rank it first among 
those studied. However, when t^ is near 1/X Procedure 3 appears to give 
comparable results with less computation. 

Procedure 3 performs very well in cases where t is approximately 
equal to l/X ; however, when t^ < < l/\ few failures tend to occur and 
the random interspersing of censored times causes the sample variance 
of the procedure to be high. In general, the empirical coverage prob- 
ability for this procedure is close to the chosen confidence coefficient. 
The actual computations needed are comparable to those for procedure 1. 
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The estimators for X defined under Procedures 2 and 3 have been 



tabulated; the one given under Procedure 2 seems to be a nearly unbiased 
estimator for X as stated. The estimator arising from Procedure 3 is 
consistently greater than the true value of X. 

The comparison of Procedures 3 and 4 are based on an abbreviated 
Monte Carlo study. This is due to the accuracy limitations of the 
computer when computing small values of X and also the added dependence 
upon the time till failure in the iterative step of the computation. 

A more extensive study will be necessary to obtain more meaningful 
results. It appears that the interval length obtained from using Pro- 
cedure 4 is not as good as that of Procedure 3; however, the empirical 
coverage probability for Procedure 4 seems to be close to the desired 
confidence coefficient and the estimator for X appears to provide a 
nearly unbiased estimator. 
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V. COMPUTATIONAL PROCEDURES 



The IBM/System 360 Model 67 computer with the Fortran IV programming 
language has been used for all computations. 

The exponentially distributed random sample is obtained by first 
generating random numbers from a uniform distribution on the interval 
(0, 1) and then by inverting these numbers to get a random sample from 
an exponential distribution with parameter X. This is done by using the 
fact that if a random variable Z has an exponential distribution with 
dis tr ibut ion function 



z < o 



-Xz 



F (Z) = P[Z < z] = J 1-e o < z < 1 

1 z > 1 

-Xz 

the rando variable X = 1-e ' has a uniform distribution on the interval 

(0, 1). Thus, if x is a random number from the interval (0, 1), the 

number z = - In (l-x)/X is a random sample number from an exponential 

distribution with parameter X. The subroutine RANDU of the IBM 

Scientific Library is used to generate these random numbers. 

The second step in the computing procedure requires a simulation of 

censored testing. Therefore, each random number z is compared with the 

pre-determined censoring time t and a counter is used to obtain the 

o 

number of failures k = . For Procedure 1, only the value of k is 

needed; however, Procedures 2, 3, and 4 require additional information. 

If the value of z exceeds the value of t^, it is disregarded and replaced 

by t ; these values of t and z which are less than t are summed to 
J o 5 o o 

obtain the observed value of the random variable XX . for use in Procedure 

l 

2. The necessity of order in Procedures 3 and 4 require the random 
numbers to be dimensioned in an array to maintain the order in which 
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they are generated. This array contains a random ordering of z f s and 

t f s which are summed to the point of the failure after the entire 

o 

sample of size n has been generated. This sum contains some of the 

truncated times t but generally not all and it contains all z < t . 

o o 

The value of this sum is the observed value of the random variable W. 

k 

of Procedures 3 and 4. 

The computation of Procedure 2 is trivial after the required normal 
deviate is read from a data card. However, due to the high number of 
possibilities for the value of k (0 < k < n) subroutines of the Beta 
and Chi-square distribution functions are used. These subroutines are 
written as cumulative distribution functions and give only the resulting 
probability, given the required input parameter. Since the input parameter 
X is unknown, it has been required to write an additional function sub- 
program to iterate toward the desired parameter given the confidence 
coefficient. The restrictions on the input parameters for these sub- 
routines force the cases of k = 0 and k = n to be ignored. This is done 
by disregarding the sample which produces k and by accounting for it in the 
value of the number of repetitions. The IBM Scientific Library subroutines 
BDTR and CDTR along with the function subprograms written by the author are 
used for these computations. 

For each of the thousand samples the confidence interval length is com- 
puted and is tested to see if it actually covers X. The final number of 
those covering X is divided by the number of repetitions yielding the 
empirical coverage probability. The sample variance is computed by using 

the relation = Qc (x-x^)^/(m-l) where x is the average length of the 

t h 

intervals, x^ is the length of the i interval and m is the number of 
repet it ions „ 



15 



Tab] e 1 CONFIDENCE COEFFICIENT 1-a » .95 





n= 






30 


40 


50 


o 


PROC. 




AVG. 


VAR. 


C.P. 


AVG. 


VAR. 


C.P. 


AVG. j 


j VAR. 


C.P. 




i. 


■ 

: 


.159 


.001 


.988 


.136 


.000 


.971 ! 


.3 20 | 


i ~ 

\ ! 

? .000 


.978 


3 


2. 




.148 


.001 ' 


.963 


.127 


.000 


.959 : 


j 

! .112 


1 .000 


.970 




3. 


i 

) 

i 


.357 


.001 


.951 


.132 


.001 


.953 


i 

S .116 ! 


| .000 


.960 




1. 


! 

1 


.107 ; 


! .ooo ! 


.961 


.090 


.000 


.977 : 


| .080 \ .000 


.963 


0 


2. 


| 




.092 | 


.ooo ; 


.951 


o 

00 

o 

• 


1 1 1 

| .000 


.951 : 


! .071 


< :ooo 


.955 




3. 


_L 


.093 ! 


• .000 ! 


.954 


.081 1 


.000 


.945 : 


( 

i .o7i : 


| .000 


.958 




1. 


] 


| 


.234 : 


: .002 * 


.978 


.200 | 


.001 


.963 j .176 


i .ooi 


.964 


3. 2 


2. 


r 

> 

t 


.213 1 


.001 ! 


.961 


.183 ! 


! .001 


.958 


i .163 


S .ooo 


.952 




3. 




.217 ; 


.ooi : 


.954 


.186 


.001 


.961 i 


l .165 ! 


> .000 


.945 




1. 


! 


.943 ! 


i j 

: .036 1 


.966 


.799 | 


.017 


.959 ; 


i .703 1 


; .011 ! .961 


82 


2. 


i 

i 

J 


! 

i 


.855 


; .023 ' 


. 942 ! 


.732 


.012 


.948 ; 


! .647 


j j 

• .007 i 


1 .945 




3. 






.871 ■ 


: .026 : 


.936 | 


.743 1 


.012 


! .954 \ .657 ' 


| .008 ! 


.936 



11 


1 . j 


! 4.54 ; .822 ! .954 j 3.84 | .423 


.966 ! 3.43 { .269 \ .957 


2. | 


4.23 ! .664 I .946 | 3.60 | .349 


.955 | 3.22 | .227 [ .954 


3. 1 


I 4.56 ! 7.10 ! .940 j 3.73 j .473 


i « j 

.949 ! 3.32 ;• .292 | .948 


33 


1. 


l * j 

3.18 ! .428 ) .970 ! 2.73 \ .227 


.972 j 2.38 1 .122 j .962 


2. 1 


! 2.77 ! .196 1 .963 j 2.41 j .113 


.956 \ 2.12 ^ .066 } .957 


3. 


2.80 i .222 \ .959 j 2.42 .118 


.957 j 2.14 ! .071 i .961 


12 


1. 


6.04 j 1.29 • .974 S 5.16 1 .704 


.962 3 4.56 l .420 1 .957 


2. 


5.50 1 .874 > .951 j 4.75 \ .499 


.955 j 4.22 t .311 ! .955 


. 3. 


5.62 \ 1.00 ’ .951 1 4.84 ? .562 


.959 1 4.28 j .338 • .955 


.06 


1 . 


12.1 f 5.39 .970 ! 10.2 \ 2.60 


.966 f 9.07 f 1.56 1 .965 


2. 


11.1 ! 3.70 - .950 | 9.39 1 1.85 


; .961 | 8.40 | 1.15 t .964 


3. 


11.4 >4.14 • .949 \ 9.55 ! 2.02 


.956 \ 8.54 l 1.29 * .962 


.1 


1 . 


10.7 ? 5.47 > .961 j 9.08 ! 2.62 | .971 1 7.98 1 1.63 i .959 


2. 


9.32 >2.46 ! .945 ! 7.99 f 1.30 ! .957 l 7.08 \ .798 \ .950 


3. . 


9.43 j 2.74 } .939 1 8.05 \ 1.32 I .958 ,* 7.12 < .827 .950 



WG. = AVERAGE LENGTH 
MR. = SAMPLE VARIANCE 
-.P. = COVERAGE PROBABILITY 
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Tabic 2 CONFIDENCE COEFFICIENT 1-cc = .99 





n= 


30 


40 


i 

| 50 






AVG. 


VAR. 


C.P. 


AVG. 


VAR. 


C.P. 


; AVG. 


VAR. 


C.P. 


t 

o 


PROC. 






















i. 


.206 


i 

| .002 


j .998 


! .177 




| .001 


1 

| .998 


.156 


; .000 


5 .996 


3 


2. 


.202 


j .001 


| .996 


.172 


| .001 


j .990 


.152 


| .000 


j .989 




3. 


.209 


.002 


.993 


.176 


j .001 


| .993 


J .153 


| .000 


.994 




1. 


.139 


.001 


.991 


.117 


j .000 


j .998 


j . 104 


1 .000 


1 .991 


10 


2. 


.123 


.000 


| .988 


.106 


j .000 


! .992 


.094 


j .000 


\ .990 

i 




3. 


.123 


.000 


| .988 


.106 


| .000 


j .995 


.094 


I .000 


) .992 


- 


1. 


.304 


.003 


! .998 


.260 


| .002 


j .996 


.232 


j .001 


j .994 


3.2 


2. 


.286 


.002 


.995 


.245 


| .001 


| .993 


.220 


! . 001 


| .986 




3. 


.287 


.002 


.994 


.246 


j .001 


| .995 


.220 


( .001 


| .989 


j 


1. 


1.22 


.060 


.988 


1.04 


I .030 

f 


J .989 


.913 


1 .018 


\ .991 


! .82 


2. 


1.15 | 


| .040 


! .987 1 

} 


| .979 


] .020 


! .988 | 

l . ! 


1 .862 


\ .013 


j .993 




3. 


1.15 


| .046 


| .988 


| .980 


j .022 


.988 


.867 


| .015 


j .988 




1. 


5.90 j 


1 1.38 


1 -996 


4.99 


| .713 


.995 


1 4.46 


j .454 


| .994 


! - 11 


2. 


5.78 | 


; l.ii 


| .989 


4.87 


| .587 


.989 | 


! 4.34 

f 


| .385 


| .990 




3. 


6.07 j 


14.6 


| .992 

i _ 


4.94 ; 


( .813 


.985 j 


■ 4.39 


j .503 


| .988 




1. 


4.13 


.729 


| .993 


3.55 ! 


| .384 


.997 i 


3.12 


| .238 


| .994 


.33 


2. 


3.70 


.342 j 


.992 


3.21 j 


| .197 : 


.994 | 


2.84 


i .132 


j .987 




3. 


3.70 


.385 


.989 


3.20 


.204 


| .994 


2.83 


| .135 


| .986 




1. 


7.84 ■ 


2.18 


.995 


6.70 


1.19 


| .995 


5.92 ' 


\ .710 

1 - 


S .997 


.12 


2. 


7.41 

1 


1.50 


.988 


6.36 


.859 


1 .991 


5.64 


.536 


j .994 




3. 


7.43 ! 


1.73 


.987 


6.39 


.970 


.992 


5.64 


.583 


j .994 




1. 


15.7 j 


9.16 


.993 


13.4 


4.68 | 


.993 


11.8 


3.05 


1 .991 

9 


.06 


2. 


14.8 | 


6.20 


.985 


12.7 


3.38 


.991 


11.2 


2.14 


! .992 

i 




3. 


14.9 j 


7.31 


.986 


12.7 


3.79 


.991 


11.2 


2.31 j 


| .993 




1. 


14.0 I 

. _ . . I 


9.33 


.989 


11.8 


4.44 


.997 


10.4 


2.49 


.998 


.1 


2. 


12.5 | 


4.28 


.981 


10.6 


2.26 


.992 


9.44 


1.34 


.996 




3. 


• 12.4 


4.74 


.983 


10.6 


2.28 


.994 


9.43 


1.44 


.994 



t 
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TABLE III 
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Avg . = Average Value of Estimate Var. = Sample Variance 



TABLE IV 



CC “ .05 n - 30 



X 


t 

O 


Proc . 


Avg. 


Var. 


C.p. 


A 

X 


Var. X 


.8 


.82 


3 


.884 


.026 


.96 


.931 


..055 


4 


1.13 


.005 


.96 


.805 


.075 


3 


. 11 


3 


4.59 


1.06 


.91 


3.92 


1.55 


4 


5.92 


4.42 


.95 


3.14 


1.30 


3 


.33 


3 


2.72 


. 192 


.95 


3.15 


.475 


4 


3.06 


.030 


.94 


2.98 


.626 


r 


.06 


3 


7.67 


3.44 


.96 


6.16 


3.98 


J 


4 


9.63 


6.41 


.98 


4.96 


3.15 


c 


.12 


3 


5.69 


1.05 


.95 


5.79 


2.07 


J 


4 


7.46 


. 123 


.95 


4.89 


2.76 


10 


.06 


3 


10.9 


3.45 


.96 


10.8 


7.44 


4 


14.4 


.642 


.98 


9.02 


9.16 



Avg. = Average Interval Length 
Vai . = Sample Variance 

c.p. = Coverage Probability 

a . 

X = Estimate of X 

A A 

Var. X = Var: nee of X 
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MAIN PROGRAM FOR COMPUTING CONFIDENCF INTERVALS, SAMPLE VARl 
ANCE, ANO EMPIRICAL PROBABILITY FOP PROCEOIIRFS 1,2. AND 3o 
NN THE NUMBER rip DIFFERENT COMBINATIONS OF CONFIDENCE 
COEFFICIENTS AND SAMPLE S I Z F 

NCOM THE NUMBER CF DIFFERENT COMBINATIONS OF LAMDA AMD 
CENSORING TIMfs 
N SAMPLE SIZE 
L° E P NUMB FR OF REPETITIONS 

IRAND INPUT PARAMETER FOR SUBROUTINE RANDU 
PROB l-ALPHA/2 ALPHA IS THE CONFIDENCE COEF® 

XNORM THE NORMAL DEVIATE 
AMD A THE PARAMETER LAMDA 
T CENSORING TIME 

DIMENSION BETA ( 101 ) , AMD! 1000 ), AMDD ( 1 COO ) , AHAT ( 1000 ) , XX 

( 100) , 

1CHI ( 100 ) , AM (1000 ) , CHIL ( 100 ) , AML ( 1000 ) , AMDI. ( 1000) , BETAL 

( 101 ) , 

1 A M DDL ( 10 00 ) , C AMD A ( 1 000 ) , CC ( 1 00 0 > , BB ( 1000) ,0( 1000) 

READ ( 5 * 1 3 ) NN 
13 FORMAT (14) 

D02 5 11 = 1 ,NN 

READ( 5, 3 )NCCM,N, LREP , I R AND , PRO B , XNORM 

3 FORMAT (417, 2X,2F10o 6) 

NRFP=LREP 

DO 1 1 0 I = 1 , N 
CHIL (2 * I ) = 0o 

110 CHI ( 2*1 ) =0o 
D0111 T = 1 , N 
RFTAL( I )=0o 

111 BETA (I ) =0 o 
PROB L= 1 o -P RGB 
D025KK=1 , NCCM 
LREP=NREP 

SC AM = Oo 
MCUN- 0 
L OUN =0 
KCUN=0 
MCUNT=0® 

AH = 0® 

L CUNT = 0 
K QUNT = 0 
CH = Co 
NC = 0 
00=0® 

N0 = 0 
B I N= Oo 
NB = 0 

READ (5 , 2 ) AMDA , T 
2 FORMAT ( 2F20o 8 ) 

D015 1 = 1 , LREP 
K=0 

SMALT=Qo 

SONOR= XNORM* XNORM 
GENERATE THE SAMPLE 
D07L =1 , N 

CALL RANOU( IRAND, I YRAN, X) 

IRAND= I YPAN 
X = - A LOG ( 1® -X) /AMDA 
IF(X-T)4, 5,5 

4 K=K+ 1 
GO TO 6 

5 X = T 

6 SMALT=SMALT+X 
X X ( L ) = X 

7 CONTINUE 
IF(K-N) 50,45,50 

45 L RE P =L R E P- 1 
1=1-1 
GO TO 15 
50 CONT INUE 

I F ( K ) 49 , 45 , 49 
49 TIMF=0o 
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oo55 I t = l » m 

IF ( XX( N- IT4 1 )-T ) 51 , 55, 55 
51 JT=N-I T+] 

GO TO 56 

55 CONT INUF 

56 00 60 I X = 1 , JT 

T IME=T IME + XX( IX) 

60 COOT IOUE 

COMPUTATION FOR PR0Co3 

I F ( C H I ( 2*K) ) 65 , 61 , 65 

61 G= 2* K 

CHI ( 2*K ) =CVAL(PROft,G) 

65 CONT IOUE 

AM ( I I = C H I ( 2 * K ) / ( 2 o * T I M E ) 

I F ( A M ( I )-AMOA) 67,66,66 

66 MCHJNT= MOUNT + 1 

67 CONTINUE 

I F(CHIL( 2*-K ) ) 165,161,165 
161 G=2*K 

CHIL(2+K)=CVAL(PR0BL,G) 

165 CONTINUE 

AML( n=CHIl(2*K)/(2n*TIME) 

IF (AML(I)-AMDA)1 66,166,167 

166 MOUN = MQUN + 1 

167 CONTINUE 

C C ( I ) = A M ( I 1 - A M L ( I ) 

CH = C H+CC ( I ) 

TF( MOUN+MOUNT-2) 171 , 170, 171 

170 NC =NC + 1 

171 NOUN =0 
MOUNT=0 

CAMOAd ) = CCMI { 2^K1 +CHIL< 2^K» >/(4.*TIME) 

SCAM=SCAM+C AMOA( I ) 

COMPUTATION FOR PR0Co2 
FN=N 

AMDAl = (K/SMALT)/(l<,+lo/(2o*FN) ) 

AH= AH+AMDA 1 
AHAT (I ) = A MDA 1 
C=U /SMALT 
B=C*C 

A MO ( I ) = ( SQNCR*C + 2o'*AMDA!+SQRT(SQN0R'*SQN0P*B + 4- o *AMDAl*S 
QNOR*C) ) / 2 

I F ( A MO ( I I - A MO A ) 9 , 8 , B 

8 KOUNT=KOUNT + 1 

9 CONT INUF 

AMQLU )=(SGN0R' t C+2o^AM0Al- 

1 SQRT( SQN0R*SQN0R*B+4 9 *AMDA1*S 

GNOR*C ) ) / 2 

IF (AMOK I) -AMD A) 108, 108, 109 

108 KOUN=KOUN+ 1 

109 CONTINUE 

0(1) =A MO ( I ) -A MDL ( I ) 

00=00+0 ( I ) 

I F ( KOUN + KOUNT-2) 173,17 2, 17 3 

172 N0=N0+1 

173 KOUN =0 
KOUNT = 0 

CCM PUT A T I ON FOR PPOCd 
K = K + 3 

IF (RET A (K) ) 20 5,2 00, 20 5 

200 B 1 =N-K+ 1 
A 1 = K 

B FT = XV A 1. ( PR OR , Al , B1 ) 

IF (BFT) 201 , 202,202 

201 GO TO 

202 CONT INUF 

BE T A ( K ) =BE T 
205 PUP= BFT A ( K ) 

A MOO ( I ) = - A LFG ( 1 o -PIJR ) / T 
IF ( A MOO ( l ) -AMO A) 11,10,10 

10 L OUN T= L OUNF + 1 

11 CONTINUE 
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I F ( B FT Al (K ) 1305, 300,305 
300 Bl=N-K+2 
A 1 = K- 1 

BET=XVAL( PRCBL, Al, Bl) 

IF(BET)?C1, 302,302 
302 CONTINUE 

BFTAKK 1=BFT 
305 PLP=BETAL(K) 

AMDOL ( I)=-ALOG(lo-PLP)/T 
I F( AMDOL ( I ) -AMDA1310, 310,311 

310 LOUN=LOUN+ 1 

311 CONTINUE 

BB ( I ) = A MOD ( I l-AMDDL ( I ) 

B IN= RI N + BB ( I ) 

I F ( L CUN + L HUNT- 2 ) 175,174, 175 

174 N B = N B + 1 

175 L GUN =0 
L CUNT = P 

15 CONTINUF 
REP=LRFP 

V AR3 =0o 
VAR1=0, 

VAR2=0, 

AVH=Oo 
SH = Oo 
00=00/ RE P 
BIN=BIN/RE D 
CH=C H/REP 
SCAM=SCAM/REP 
AH=AH/LRFP 
0020 J=1 , LREP 

SH = SHHC AMO A( J ) -SCAM )»72 
AVH=AVH+ ( AM AT ( J ) -AH)**2 

V A R 1 = V A R 1 ( 0 ( JJ-OD)-^? 

VAR2=VAR?+( PB(J)-BIN)^*2 
VAR 3 = VAR 3+ ( CC (JI-CH}**2 

20 CONTINUE 
PB=NB/REP 
PC=NC/RFP 
PO=NO/REP 
SH=SH/ ( RFP-lo ) 

VAR1=VAR1/ ( R F P-1 0 ) 

VAR2=VAR2/ (REP-1 ) 

VAR3 = V AR3/ ( REP-1, ) 

AVM=AVH/ (REP-lo ) 

WRITFI6, 21 ) 

21 FORMAT ( IX, 'REPS' ,3X, ‘SAMPLE SIZE LAMDA CUTOFF TIME PRO 
BAB IL ITY' ) 

W R I T E ( 6 , 2 2 ) LRFP,N, AMOA, T, PROB 

22 FORMAT ( // , 2 I 8 , 1 X, 3 F 10o 5 , // I 
WRITF(6,?3) 

23 FORMAT ( 14X , »0o Do P ROC ' , 7 X , » B I N „ MF TH • , 6 X , • CH I - SO, METH 

' ) 

WRITF( 5, 26 ) 00 ,BI N,CH, VAR1 , VAR2 , VAR 3 , PO, PB, PC 
26 FORMAT ( IX, ' AVERAGE' ,3( 5X,F10o5 ),//, 1 X, 'VARIANCE* ,3(5X, 
F 1 Do 5 ) , //, 

1 1 X , 'COVERAGE PROB' , 3 ( 5 X , FI 0, 5 ) ) 

WRITF(6 ,24 ) AM, AVH 
WR I TE ( 6 , 24 ) SC AM f SH 

24 FORMAT ( IX, ' AVG, LAMDA HAT=' ,F12« 6, 'VARIANCE LAMDA H AT = 
• » F 1 2 , 6 ) 

25 CONTINUF 
FND 

XVAL ITERATES TOWARD LAMDA FOR PROCEDURE 1 



FUNCTION XVAL( ALPHA, Al ,B1) 
DATA FPS/o 5F-07/ ,N/1000/ 

A = A 1 
R = B 1 

A L = A L P H A 
4 XLAST=- 0 1 
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M = 0 
X = o 1 

0 X =« 1 

11 CALL BDTP( X,A,D,P, D, IFR) 

TEST=P-AL 

TM TEST) 13,5,13 

13 IF( ARS( TFST >- o 5E-04) 12,12,6 

12 IF(A3S( ( T EST-XLAST)/XLAST)-EPS ) 5 , 5 , 6 

5 XVAL = X 
RETURN 

6 I F< TEST*XL AST ) 7, 8, 8 

7 DX=-„5*DX 

8 X=X+DX 
XLAST=TEST 
M = M+ l 

1 F ( M -N ) 1 1 , 9 , 9 

9 WRITE (6,15) 

15 FORMAT( IX, • NO CONVERGENCE') 

XVAL=-lo 

RETURN 

END 

CVAL ITERATES TOWARD LAMDA FOR PROCEDURE 3 



FUNCTION CVAL ( ALPHA, G) 

DATA EPS/. 5E-G7/ ,N/1000/ 

AL=ALPH A 
A XLAST=-ol 
M=0 
X = o 1 
DX = . 1 

11 CALL COT R ( X , G , P, D, I ER ) 

TE ST = P-AL 

I F ( TEST ) 13,5,13 

13 IF(ABS( TEST )-.5F-04) 12,12,6 

12 1F(ABS( ( TEST- XL A ST) /X L AST ) -E PS ) 5 , 5 , 6 

5 C V A L - X 
RETURN 

6 IF(TEST*XLAST)7,8, 8 

7 DX=-o5 J OX 

8 X=X+DX 
XLAST=TFST 
M= M+ 1 

I F ( M-N ) 1 1 , 9 , 9 
o WRITE! 6, 15) 

15 FORMAT ( IX, 'MO CONVERGENCE') 

CVAL=- lo 

RETURN 

END 
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